/******************************************************************************

	COPYRIGHT(C) JONAS 'SORTIE' TERMANSEN 2009, 2010.

	This file is part of the Maxsi Library.

	Maxsi Library is free software: you can redistribute it and/or modify it
	under the terms of the GNU Lesser General Public License as published by
	the Free Software Foundation, either version 3 of the License, or (at your
	option) any later version.

	Maxsi Library is distributed in the hope that it will be useful, but
	WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
	or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
	License for more details.

	You should have received a copy of the GNU Lesser General Public License
	along with Maxsi Library. If not, see <http://www.gnu.org/licenses/>.

	Maxsi Library
	A powerful Cross-Platform C++ General Purpose Library that allows you to
	efficiently create high-performance and reliable applications.

	MaxsiRegression.cpp
	Regressions using the least squares method. Mathematically proven! And
	Patent-Free! And Abestos-Free! And Software-Free!

******************************************************************************/

#include "MaxsiLibrary.h"
#include <math.h>

BeginMaxsiNamespace

//=============================================================================
// y = ax + b
// For a dataset with n values of X and Y, we get that
// a = (sum(i=1, n, Xi * Yi) - n*avg(X)*avg(Y)) / (sum(i=1, n, Xi^2) + sum(i=1, n, Xi*avg(X))
// b = avg(Y) - a*avg(X)
// I have a wonderful proof for this, but there is not room for it here. <Sortie@Maxsi.dk>
//=============================================================================
bool LinearRegression(size_t N, double* X, double* Y, double* A, double* B)
{
	if ( N == 0 || X ==  NULL || Y == NULL || A == NULL || B == NULL ) { return false; }

	double	SumXiYi		=		0.0;
	double	SumXiXi		=		0.0;
	double	AverageX	=		0.0;
	double	AverageY	=		0.0;

	for (size_t I = 0; I < N; I++)
	{
		SumXiYi			+=		(X[I])*(Y[I]);
		SumXiXi			+=		(X[I])*(X[I]);
		AverageX		+=		X[I];
		AverageY		+=		Y[I];
	}

	AverageX			/=		N;
	AverageY			/=		N;

	*A					=		(SumXiYi-N*AverageX*AverageY)/(SumXiXi-N*AverageX*AverageX);
	*B					=		AverageY - (*A) * AverageX;

	return	true;
}

//=============================================================================
//	Same as above, but with b = 0.
//=============================================================================
bool ProportionalRegression(size_t N, double* X, double* Y, double* A)
{
	double	AverageX	=		0.0;
	double	AverageY	=		0.0;

	for (size_t I = 0; I < N; I++)
	{
		AverageX		+=		X[I];
		AverageY		+=		Y[I];
	}

	//AverageX			/=		N;
	//AverageY			/=		N;

	*A					=		AverageY / AverageX;

	return	true;
}

//=============================================================================
//	I don't have a proof for this but it seems to generate the same values as
//	other programs, so we're all good!
//=============================================================================
double Correlation(size_t N, double* X, double* Y)
{
	double	AverageX	=		0.0;
	double	AverageY	=		0.0;

	for (size_t I = 0; I < N; I++)
	{
		AverageX		+=		X[I];
		AverageY		+=		Y[I];
	}

	AverageX			/=		N;
	AverageY			/=		N;

	double	Upper		=		0.0;
	double	Lower1		=		0.0;
	double	Lower2		=		0.0;

	for (size_t I = 0; I < N; I++)
	{
		Upper			+=		(X[I]-AverageX)*(Y[I]-AverageY);
		Lower1			+=		(X[I]-AverageX)*(X[I]-AverageX);
		Lower2			+=		(Y[I]-AverageY)*(Y[I]-AverageY);
	}

	Lower1				=		sqrt(Lower1);
	Lower2				=		sqrt(Lower2);

	return	Upper/(Lower1*Lower2);
}

EndMaxsiNamespace

